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We investigate by numerical simulations and analytical calculations the Bak-Sneppen model for 
biological evolution in scale-free networks. By using large scale numerical simulations, we study the 
avalanche size distribution and the activity time behavior at nodes with different connectivities. We 
argue the absence of a critical barrier and its associated critical behavior for infinite size systems. 
These findings are supported by a single site mean-field analytic treatment of the model. 
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'"q . Many real systems, ranging from biological systems such as food webs |1 3fl to communication systems [|4|,|5|], exhibit 
properties that lies in between those of regular lattices and random graphs pf. They are usually referred to as complex 
networks . These networks may have large clustering coefficients like regular lattices but also have a small diameter 
which is a typical feature of random graphs. Among the class of complex networks, a particular role is played by 
scale-free networks (SF) || in which there are not characteristic fluctuations in the connectivity of the nodes. This 
implies that the probability P^ that a node is connected with other k nodes follows a power law, i.e., P^ ~ fc -7 , in 
! t \ contrast to exponential graphs in which P^ is exponentially bounded [^0^1 . 

SF networks have been recently recognized to describe several real growing networks p|-^ Jj|jL0[ , and at the same 
time have proved to show very peculiar features with respect to physical properties such as damage tolerance [p"T|| , 
epidemic spreading [Jl3, and diffusion properties p4[ . It is then natural to ask whether and to what extent the 
topology of these complex networks would affect many of the results obtained for punctuated evolution models in 
regular lattices. 

In this Letter, we study the Bak-Sneppen model (BS) jl5| on SF networks. We perform large scale numerical 
simulations and find that contrary to what is observed in regular lattices and in exponential networks, the system self- 
organizes into a stationary state characterized by the lack of a critical threshold barrier in the thermodynamic limit. 
This result is confirmed analytically by constructing the geanological tree of the avalanches and performing a mean- 
field approach which takes into account the strong fluctuations of the network's connectivity distribution. Finally, we 
discuss the consequences that the present results could have in the general context of punctuated evolution models. 

The standard BS model psj can be considered as an ecological system formed by many species that interact one 
with each other if they are in contact, i.e., the interactions are local. To each node of the graph we allocate a random 
fitness barrier fcj € [0, 1) (i = 1, . . . , N) that represents the ability of species to survive or mutate The fitness 

are initially uniformly distributed between and 1 and the dynamics is updated according to the following rules: i) 
at each time step, the species with the minimum barrier b m i n is located and mutated by assigning a new random 
value for its fitness. This somehow mimics the Darwinian principle that the least fit species evolve by mutation or 
disappear; and ii) the species directly linked to the species with b min change their fitnesses to new random numbers as 
the result of their interactions. By applying these rules repeatedly, after a transient period and regardless of the initial 
distribution of fitness barriers, the fitness distribution P(b) of the system evolve toward a self-organized stationary 
state where several physical quantities can be measured. In previous studies on regular geometries and exponential 
networks |l6|,[l^,[l9|] , it has been found that the distribution of barriers tends, in the limit of infinite system sizes, into 
a step-like function characterized by the existence of a single parameter b c . In these cases, for b < b c the distribution 
P(b) = 0, otherwise, P(b) is equal to a constant value. 

In order to study the BS model on a SF topology we use the network obtained using the algorithm of Ref. || . This 
is a stochastic growing network model in which at each time step a new node (or vertex) is added to the network 
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and connected preferentially to the already existing ones with a probability that depends on their connectivities. In 
practice, we start with a small number mo of disconnected nodes and at each time step the network grows by adding 
a new vertex. This is connected to m old ones i with a probability Tl(ki) — kj / Yl . kj. By iterating this scheme, 

a network of size N with connectivity distribution ~ fc~ 3 develops. As a first step, we performed numerical 
simulations of the BS model in SF networks with sizes ranging from N = 10 3 to N = 5 x 10 5 . The BS dynamics 
is iterated to achieve a stationary state giving rise to the step function behavior usually observed in regular lattice. 
On SF networks, however, it appears that the critical barrier b c is not an intrinsic quantity and that b c — > if 
N — > oo. Figure [l] shows the dependency of the critical barrier with the system size. It turns out that as the system 
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size is increased the critical barrier b c shifts leftward and goes to zero logarithmically as b c ~ 1/lnN. Thus, in the 
thermodynamic limit (N — > oo), there is no threshold barrier; i.e. b c — 0. Above the threshold, the distribution of 
barriers is flat with overimposed statistical fluctuations. A further check of this behavior is provided by the analysis 
of the burst activity (avalanche) in the system. An avalanche is defined as the number of subsequent mutations below 
a certain threshold. In Fig. |^, we show the avalanche size distribution for two different values of avalanche threshold 
b. For this system, which consists of N — 10 5 species, the critical barrier is b c — 0.044. For values of b below b c , 
the avalanches are distributed according to a power law P(s) ~ s~ T with r = 1.55 ± 0.05, a value close to the mean 
field exponent r = 3/2 observed in other versions of the BS model jl7| and in self-organized critical models Jl8|]. The 
inset illustrates that as soon as the 6-barrier is placed above the threshold, the distribution splits into two parts (see 
the inset) and is characterized by the excess of large avalanches signaling that we are in the supercritical region. If 
the barrier 6 were further increased, then there would be only one avalanche whose size will be limited only by the 
observation time. Besides, for very large system sizes the absence of a critical barrier makes that even for small b the 
avalanche size distribution looks like the inset of Fig ^ indicating the lack of a critical point. 

The numerical evidence tells us that the BS model behavior radically changes in SF networks. The lack of a 
threshold can be intuitively understood by recalling that for regular lattices and for networks in which the connectivity 
distribution is exponentially bounded, the threshold barrier decreases as the number of neighbors with which a species 
interacts increases. This is the case studied in |Ts|| , where an ecology consisting of N species, each one interacting 
with Zi neighbors, was studied. In this model, the z% neighbors are drawn from a Poisson distribution with a mean (z). 
When increasing (z) and thus the connectivity, it was observed that the threshold barrier b c decreases. In SF networks 
like the one studied here, the fact that the fluctuations in the number of neighbors of each species diverge ((k 2 ) = oo) 
makes null the threshold barrier when N — > oo. It is worth remarking that the absence of a critical threshold due to 
the strong fluctuations of the SF networks connectivity has been recently reported for epidemic spreading [ p^|JT^ ] and 
percolation- like phase transitions . 

In order to support with analytical considerations the absence of a critical barrier, we will analyze the burst like 
activity observed in the evolution process following a MF approach. MF approaches are expected to give the right 
solution for growing random networks which are defined by non-local random topologies. Let us define a b avalanche 
when there is one specie i with bi < b while for all the others bj > b. The b avalanche lasts until all the barriers are 
above the b barrier. Thus, if there is a critical barrier & c , the avalanches should always be finite for b < b c since below 
the threshold it is expected that no species remain forever. Otherwise, there is a nonzero probability to observe an 
infinite avalanche that approaches one as we move away from the critical barrier. Let us assume that < b < b c . 
Thus, a b avalanche must be finite. 

Within a b avalanche one can build a geanological tree as it was shown in Refs jnj . A node in this tree represents 
a species whose barrier goes below b at some step. Since we have assumed that the avalanche is finite this species will 
be selected for evolution at some later step, assigning a new random barrier to it and its neighbors. In this process, 
the barrier of the species itself and those of its neighbors may go below the threshold again and, in such a case, they 
are represented by other nodes in the tree. The fact that these nodes are causally related is represented by a direct 
link from the ancestor (the species with the minimum barrier) to its sons (the species whose barriers go below 6). 
The generation t of a node in this tree is then defined as the number of nodes one needs to pass in order to arrive to 
that node starting from the node that triggered the avalanche. As we have assumed that the avalanche is finite, the 
relative density of species whose fitness values are below the critical barrier should vanish in the stationary state. 

Let <j>k{t) be the relative density of nodes at generation t with given connectivity k. This connectivity is taken as 
the connectivity of the real SF network, where the BS evolution rules take place, and not to that corresponding to 
the geanological tree. The rate equations for fa (t) are given by 



The first term on the r.h.s. takes into account that, when a species is selected for evolution, with probability 1 — 6 
its new randomly assigned barrier may be above b. The second term considers the fraction of nodes 1 — fa(t) whose 
barriers are above b, but change their barriers because they have a link to the species selected for evolution. The 
factor b takes into account that the new barrier will be below b with probability b. 0(b) denotes the probability that 
a link points to a node selected by evolution, and k@(b) the probability that a species with connectivity k have a 
link pointing to the species which are going to evolve. By using the MF considerations of Ref. p2[ , this probability 
is given by the number of links belonging to species with connectivity fc, kPj-N, divided by the total number of links 
'^2 S sP s N. Hence, making the summation over k it results that 



d t fa(t) 



(l-b)fa(t)+bk[l-fa(t)]Q(b). 
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In the stationary state dt4>k(t) = and from Eq. (|y) and (|^) we obtain two self-consistent equations from which 
one obtains 0(b) and 4>k- Finally, the stationary density of nodes with barrier below b is = J^k This set of 
equations is analogous to the one studied in [fl2f . The solution gives that for SF networks with 2 < 7 < 3 and for 
any value of b the stationary density of species below b is finite and, therefore, in the thermodynamic limit there 
is a finite probability to obtain an infinite avalanche where strictly speaking, (f> is nonzero. However, this result is 
in contradiction with our initial assumption that < b < b c , for which the avalanches should be finite. Hence, we 
conclude that for 2 < 7 < 3 there is not finite threshold b c when N — » 00. For any finite system size, the threshold, 
although very small, is not zero since there is a finite probability that all species have their fitness values above the 
critical barrier at the same time. It is worth noting at this point that the above arguments can be used to easily 
show that for exponential networks with an average connectivity (k) the scenario at the stationary state is completely 
different. In this case 6(6) ~ (f> and it results that there is a finite threshold barrier which is given by b c = ppm , as 
can be seen by direct comparison with the numerical values reported in p9[ , b c = 0.3446 « 1/3 and b c = 0.2575 w 1/4 
for (k) = 2 and (k) = 3, respectively. Thus, for exponential networks, the threshold barrier is in general finite and 
it is determined by the first moment of the degree distribution (k). On the contrary, in SF networks, the lack of a 
threshold barrier takes place when the second moment (& 2 ) diverges, i.e. when 7 < 3. The peculiar topology of these 
networks gives rise to the existence of a few nodes with a very high number of neighbors, and thus the fluctuations 
in the number of neighbors emanating from each node are unbounded (that is, the second moment diverges). This 
fact becomes determinant for the dynamics of the system. We would like to remark that our results are on the same 
line as those obtained in spreading dynamics on SF networks where the critical parameter is related to the first and 
the second moments through the relation (k)/(k 2 ) [ p^JT^ ] which implies a vanishing threshold in the thermodynamic 
limit. However, we emphasize that real networks have always a finite size N and thus a finite effective threshold, 
although extremely small for very large system sizes. 

The peculiar features of the BS model depend on the highly heterogeneous nature of SF networks. As a further 
evidence of this heterogeneity, we explore how the activity in the network is distributed according to the nodes' 
connectivity. We have measured the number of times rik a node with connectivity k has been selected to mutate because 
it has the minimum barrier among all the species. The results obtained shows that the activity is concentrated in the 
species with high connectivities. This indicates that the activity patterns are strongly correlated with the connectivity 
of each site. Species with high degrees of interaction with their neighborhoods are thus more susceptible to mutate 
or become extinct than their counterparts. Speciation mechanisms p(| during evolution have been suggested as key 
natural processes to explain the evolutionary patterns found in all lineages. For speciation to be accomplished, it has 
been proposed both geographical isolation and environmental stress as important triggering factors. In the model 
studied here, the nodes in the network represent species, whereas the links between them correspond to the complex 
dependency relationships that can be established in an ecosystem, for example, through a food web or predator-prey 
relationship. Therefore, the species with high connectivities may be thought of as being more stressed than others 
with low connectivities are. Consequently, we expect that the punctuated equilibrium pattern be manifested in a 
significant reduction in the cumulative activity of the lowly-connected species as compared to that of the highly- 
connected ones. Figure || shows the cumulative activity for species with connectivity k — 50 and k = 3 (inset). As 
can be clearly seen, the punctuated equilibrium behavior appears in both sub-populations. However, the higher the 
species connectivities are, the shorter the periods of intensive activity are. Thus, for the same time window, speciation 
and rapid episodes of major evolutionary changes are expected to take place in the more stressed, highly dependent 
species. This correlation has been observed in several studies |2l| . Besides, the fact that speciation can occur rapidly 
over a large area containing millions of individuals has been documented by Williamson p2|] , who studied a series of 
mollusk fauna of the eastern Turkana basin, in Africa. 

In summary, we have studied the BS model on scale-free networks. We have shown that the highly interacting 
species play a key role in the evolution process. Numerical results as well as analytical arguments point out the 
absence of a critical threshold barrier in the thermodynamic limit. For real system, which are always of finite size, 
the threshold, although very small, exists. We have also found that the activity patterns are strongly correlated with 
the topology of the network. Finally, we have outlined the possible implications our results may have for speciation 
events. 

We are grateful to A. Vespignani for helpful comments and discussions. 
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FIG. 1. Threshold barriers for several system sizes ranging from N = 10 3 to N = 5 x 10 5 . The inset shows the distribution 
of barriers in the stationary state for systems of size 10 4 (circles) and 10 5 (squares) for b < 0.3. 
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FIG. 2. Avalanche size distributions for 6 = 0.04 and b = 0.09 (inset). The system, which consists of N = 10 5 species, has 
the threshold value at b c = 0.044. The data were recorded after a transient of 10 6 mutations. A full line corresponding to a 
power law with exponent r = 1.5 has been drawn for comparison. 
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FIG. 3. Cumulative activity for species with k = 50 and k — 3 (inset) neighbors. In this case, k = 3 corresponds to the 
minimally connected species. The activity strongly depends on the species' connectivity and for highly connected species it is 
~ 100 times that of the minimally connected sites. The punctuated equilibrium behavior can be nevertheless observed for both 
sub-populations. 
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